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Introduction 

Many studies have been carried out during the last decade to study the effect of cli- 
mate change on crop yields and other key crop characteristics. In these studies, one or 
several crop models were used to simulate crop growth and development for differ- 
ent climate scenarios that correspond to different projections of atmospheric CO 2 
concentration, temperature, and rainfall changes (Semenov et al., 1996; Tubiello 
and Ewert, 2002; White et al., 201 1). The Agricultural Model Intercomparison and 
Improvement Project (AgMIP; Rosenzweig et al., 2013) builds on these studies with 
the goal of using an ensemble of multiple crop models in order to assess effects of 
climate change scenarios for several crops in contrasting environments. 

These studies generate large datasets, including thousands of simulated crop 
yield data. They include series of yield values obtained by combining several crop 
models with different climate scenarios that are defined by several climatic variables 
(temperature, CO 2 , rainfall, etc.). Such datasets potentially provide useful informa- 
tion on the possible effects of different climate change scenarios on crop yields. 
However, it is sometimes difficult to analyze these datasets and to summarize them 
in a useful way due to their structural complexity; simulated yield data can differ 
among contrasting c lima te scenarios, sites, and crop models. Another issue is that 
it is not straightforward to extrapolate the results obtained for the scenarios to alter- 
native climate change scenarios not initially included in the simulation protocols. 
Additional d 3 maniic crop model simulations for new climate change scenarios are 
an option but this approach is costly, especially when a large number of crop models 
are used to generate the simulated data, as in AgMIP. 

Statistical models have been used to analyze responses of measured yield data 
to climate variables in past studies (Lobell et al., 201 1), but the use of a statistical 
model to analyze yields simulated by complex process-based crop models is a rather 
new idea. We demonstrate herewith that statistical methods can play an important 
role in analyzing simulated yield data sets obtained from the ensembles of process- 
based crop models. Formal statistical analysis is helpful to estimate the effects of 
different climatic variables on yield, and to describe the between-model variability 
of these effects. , 

These statistical methods can also be used to develop meta-models, i.e., statistical 
models that summarize process-based crop models, enabling scientists to explore 
the effects of new climate change scenarios. This approach is illustrated with two 
simulated yield datasets obtained by AgMIP for maize and wheat diat were generated 
by using ensembles of process-based crop models. The }deld datasets were used 
to develop a meta-model that provides a simplified representation of the original 
ensemble of crop models. The proposed meta-model is a statistical regression with 
random coefficients that describe the variability of the simulated yield data across the 
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original crop models. Once the statistical regression has been fitted to the simulated 
5 deld datasets, the meta-model can predict the effects of climate changes among the 
crop models and can thus be used to study alternative climate change scenarios fiiat 
were not initially simulated. Finally, the results obtained with the meta-model were 
used to assess the risk of yield loss for maize and wheat as functions of temperature 
and CO 2 concentration change. 

Materials and Methods 
Simulated yield data 

Yield data were simulated with 19 maize crop models and 26 wheat crop models in 
four contrasting sites for each species located in France (Lusignan), USA (Ames), 
Brazil (Rio Verde), and Tanzania (Morogoro) for maize, and in the Netherlands 
(Wageningen), Argentina (Balcarce), India (New Delhi), and Australia (Wongan 
Hills) for wheat. 

Several climate scenarios were considered. For both species, four temperature 
changes (-fO. -1-3, -|- 6 , -|-9®C) and five atmospheric CO 2 concentration changes (+0. 
-t-90, -1-180, -i-270, -1-360 ppm; where the baseline CO 2 concentration was set at 360 
ppm) were considered over the 1980-2010 time-period. In addition, two precipita- 
tion changes (0% and —30%) were considered for maize and three rates of nitrogen 
fertilization (50%, 1(X)%, and 150% of reference rates) were considered for wheat 
The other scenarios were defined by combining temperature, CO 2 concentration, 
and precipitation changes. Thirty years of yield data were generated with each crop 
model for each of the climate scenarios (Figs. 1 and 2) and the simulated yield values 
were averaged over the years. Details of the maize and wheat protocols can be found 
in Bassu et al. (2014) and Asseng et al. (2013), respectively. 

Statistical model 

Simulated maize and wheat yield data were analyzed using two-level statistical 
random-effect models (Rnheiro and Bates, 2000), Hie following statistical model 
was used to analyze the maize yield data for each site separately: 

Level 1, within crop model 

Yij = aoi + CKi/ATy -f a2iAT^ + asiACy 

-|-Q:4iACy -I- asiACijATij -h aeiARy -|- Sy (1) 

where Y/j is the 30-year mean yield simulated for the ith crop model (where i = 
i i^andfor the 7 th scenario (where 7 = 1, . .. , QO-Thevariables A7y, AC,-, 
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Fig, 1, Maize yields simulated by 15 crop models in Lusignan, France for different levels of tem- 
perature change (each letter corresponds to one crop model). The reported yield data were obtained 
for a change of CO 2 concentration equal to +180 ppm and for the baseline (1981-2010) level of 
precipitation. Only 15 out of the 19 available crop models simulated this level of [GO 2 ] and were thus 
included in the figure. 

and A/?y are the respective changes from the 1980-2010 baseline of temperature, 
atmospheric CO 2 concentration, and precipitation for model i and scenario j. The 
value of ey is assumed to be independently and identically distributed such that 
gy ~ iV(0, T^) where is the variance describingthe residual error. The coefficients, 

au (where k = 0, . . . , 6), are seven random regression coefficients distributed 
according to independent Gaussian probability distributions. 

Level 2, between crop models 


Uki ~ Nipk, of), k = 0 6 (2) 

The distribution parameters, /u.jt (where ^ = 0, . . . , 6),arethesevenmeanregression 
coefficient values that represent the mean yield baseline (/^o) and mean effects of 
temperature, CO 2 , temperature-C02 interaction, and rainfall (jui, . . . , fi^) over the 
P crop models. The seven variances, <Tj (where A: = 0, . . . , 6), describe the between- 
model variability of the random regression coefficients. 
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Fig. 2. Wheat yields simulated by the 26 crop models in Wageningen, the Netherlands, for different 
levels of temperature change {each letter corresponds to one crop model). The reported yield data were 
obtained for a change of CO 2 concentration equal to +180 ppm and for the baseline (1981-2010) 
level of fertilization. 

This statistical model assumes that the ensemble of P crop models is a gamplp. 
taken within a population that includes all possible crop models for a given crop while 
flexibly allowing for the incorporation of additional crop models in the future. The 
probability distributions defined by Equation (2) describe the between-crop-model 
variability of the yield response to climate change within the whole population of 
crop models. These probability distributions cover the ranges of climate effects 
considered in different crop models. The relationship defined by Equation (1) is 
assumed to be valid for all crop models, but its parameters, am, are assumed to vary 
among crop models. However, this statistical model only describes 30-year mean 
yield responses and is not intended to describe the year-to-year variability of crop 
yields. Creating a statistical model that captured year-to-year variability would have 
required extra random terms and additional parameters, which would have overly 
complicated the calculated model and thus was not considered here. 

The same form of model was used for wheat as was used for maize, but replaced 
the rainfall variable, AJtij, with a variable tiiat describes the level of N fer tilisati on 
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in Equation (1). For the wheat pilot studies, precipitation was held constant in the 
wheat simulation protocol. 

The population parameters of the statistical model fik, of, and were estimated 

by restricted maximum likelihood. The model-specific regression coefficients, 
(where k = 0, ... ,6 and i = 1, ..., P), were estimated by using the best linear 
unbiased predictor using the R software package “nlme” (Pinheiro and Bates, 2000) 
and the estimated values will be henceforth referred to as agiu. The model was fitted 
for each crop and each site separately, but was calculated for all crop models together. 
Results were analyzed site by site. 

Assessment of the statistical model 

The statistical model was compared to other statistical models, including models 
with fewer explanatory variables, models with fewer random coefficients, and a 
model that includes no random coefficients (i.e., classical linear regression). All the 
models were compared by using the Akaike information criterion (AIC; Pinheiro 
and Bates, 2000), where a lower AIC value corresponds with the better model, and 
we found that the model defined by Equations (1) and (2) led to lower AICs than 
the simpler models. The assumption that the residual errors, eij, were independent 
was assessed by developing another statistical model that incorporated correlated 
residual errors. The AIC value of this model was higher and the estimated cor- 
relation coefficients were very low (7.6 x lO"'* and —4.7 x 10“^ for maize and 
wheat, respectively). The statistical model’s quality of fit was also assessed by using 
graphical analysis. Evaluating the residuals did not reveal any obvious biases due 
to the assumptions made on the relationships between yields and input factors. 

Estimation of the effect of climate change on yield 

The statistical model described above was used to compute three different types of 
outputs: 

• The average yield loss/gain due to climate change over the ensemble of crop 
models. 

• The yield gain/loss estimated for individual crop models due to changes in climate 
variables. 

• The probability that the yield loss/gain will exceed a given threshold. 

For maize, the average yield difference obtained between a given climate change 
scenario (characterized by AT, AC, and AR) and the baseline scenario was 
expressed as: 

AT = ii\A.T + ^sAC -|- -f- -1- jUgAR (3) 
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The yield difference described in Equation (3) is averaged over all crop models- 
this difference corresponds to an average yield gain or to an average yield loss over 
the 19 submitted crop models. Equation (3) defines a meta-model that simnlgtes the 
average output of the original ensemble of crop models. This meta-model enables 
the computation of the yield differences for any changes in temperature, CO 2 , and 
precipitation, AT, AC, and AR, at each of the four considered sites. 

For a given crop model i, die yield difference was expressed as: 

AYi= UeuAT + tteiiAT'^ -f Q!e3/AC -f ae4,AC^ -I- UeSiACAT + a^f^AR (4) 

Equation (4) defines a meta-model simulating the output of the ith crop model. 
The yield difference (4) is crop model-specific. It corresponds to the climate ebnng i- 
effect on 3 deld that would have been obtained vrith the ith crop model if this crop 
model were run for a climate change scenario that were characterized by AT, AC, 
and AR. 

The statistical model defined by Equations (1) and (2) was also used to com- 
pute the probability that the yield difference will exceed a threshold and, more 
specifically, the probability of yield gain or of yield loss that results from a change 
in the temperature, CO 2 concentration, and/or precipitation. This probability was 
computed from the following Gaussian probability distribution N{(Iay, o\y) 

Hay = H\AT + hi^T'^ + + ii^AC^ + fisACAT + fieAR (5) 

oIy = ct\AT'^ -I- alAT^ + ofAC^ + o^AC^ -j- of AC^AT^ -|- of A/?^ (6) 

Similar relationships were derived for wheat by replacing precipitation change 
with fertilization change in Equations (3) to (6). Note that the variance defined 
by Equation (6) is not constant but varies as a function of the climate scenario 
characteristics. 

Results and Discussion 
Increases in temperature 

Figure 3 shows the statistical meta-models’ quality of fit calculated for maize and 
wheat (one site per crop). The distributions of the residuals for both Tn^t7P. and 
wheat are symmetric and do not display trends. The fitted values obtained with 
the statistical meta-models are in most cases close to the )delds Rimn1at<»d by the 
processed-based crop models. The residual standard error of the fitted statistical 
models is equal to 0.19 and 0.22 t/ha for maize and wheat, respectively, in the two 
sites considered in Fig. 3. A similar quality of fit was obtained for the other sites. 

Figure 4 shows the change in the mean yield from the baseline for one maize site 
(Fig. 4a) and one wheat site (Fig. 4b) as affected by a CO 2 concentration increase 
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(a) Maize in France (Lusignan) 
(residual standard error = 0.19 t/ha) 
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(b) Wheat in the Netherlands (Wageningen) 
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Fig. 3. Residuals of the statistical model and simulated versus fitted values obtained for (a) maize, 
and (b) wheat (one site per crop, all scenarios). The graphics compare yield values simulated by the 
crop models and yield values fitted by the statistical model for all scenarios. 


of 180 ppm ([CO 2 ] = 540 ppm) and an increase in mean seasonal temperature 
ranging from O^C to 6^’C. Each emulated model yield as calculated by using the 

crop model-specific coefficients (where k = 0 6 and i = 1, . . . , F) is 

plotted with a blac^ line and thus can be seen as a proxy of a given crop model 
Negative yield differences can be interpreted as mean yield losses and, conversely, 
positive yield differences can be interpreted as mean yield gains without the need 
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(a) Maize in France (Lusignan) 
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(b) Wheat in the Netherlands (Wageningen) 
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for re-running the original, process-based crop model. The continuous red curve 
indicates the mean yield response to the given climate scenario as compared to 
the baseline, i.e., the effect averaged over all crop models. The red dashed curves 
indicate the first (probability = 0.25) and third (probability = 0.75) quartiles of the 
climate change effect. 

According to Fig. 4, most of the crop models estimate that a temperature increase 
negatively impacts yield projections of both maize and wheat in all locations. But 
this effect is highly variable among crop models, with some models predicting little 
response to temperature. For maize. Fig. 4a illustrates how, on average across the 
ensemble of crop models, the statistical meta-model emulates a yield loss when tem- 
peratures exceed -I- 1 .5°C with a CO 2 concentration increase of -H 1 80 ppm in Lusig- 
nan, France. Contrastingly, the models suggest that wheat is more resilient to heat 
stress and would require a temperature threshold equal to -1-3.6®C before experienc- 
ing yield losses below the baseline in Wageningen, the Netherlands (Fig. 4b). Based 
on the emulated projections seen in Fig. 4b, the models indicate that, should climate 
change be limit ed to a -|-3°C change and include a CO 2 concentration increase of 
-f-180 ppm, moderate temperature increases could, in fact, lead to gains in wheat 
production in this location. 

Figure 4 also demonstrates the large variability between crop models and dis- 
plays how this variability increases as a function of temperature. The differences 
between the first and third quartiles are much larger for higher temperature changes, 
at the given CO 2 concentration. This result indicates that the differences between 
crop models are much larger for high-temperature changes than for low-temperature 
changes. The level of divergence between the predictions of the crop models thus 
tends to increase with the temperature value considered in the climate change 
scenarios. 


Increases in atmospheric CO2 concentration 

Figure 5 shows the effect of climate change at one site for both maize and wheat 
yields under increasing levels of CO 2 concentration, ranging from 0 to -1-360 ppm 
from the simulated baseline concentration (i.e., 360 to 720 ppm) at a constant temper- 
ature increase of -|-2°C. Figure 5a illustrates how all crop model simulations predict 
a positive effect on maize yields in Lusignan, France when exposed to increased 
CO 2 concentrations, but this effect is small and varies among crop models. The 
mean curve suggests that the effect of a temperature increase of -|-2°C is negative 
for the range of CO 2 concentration increases, thus the benefit of the increased CO 2 
concentration does not outweigh the negative effect of the temperature increase. On 
average, the crop models predict a yield loss for maize, even at high CO 2 concen- 
trations, due to this -f2°C temperature change. 
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(a) Maize in France (Lusignan) 



AC02 (ppm) 


(b) Wheat in the Netherlands (Wageningenj 



AC02 (ppm) 


Fig. 5. The effect of climate change on yield for a temperature increase of +2®C, and for levels of 
CO2 concentration increase ranging from 0 to +350 ppm. Results are reported for (a) and (bl 
wheat for one site per crop. Yield difference = yield obtained with a climate change scenario — yield 
obtained with the baseline (1981-2010). Each black curve corresponds to a given crop model The 
black curves were obtained widi the statistical model by using (he crop model-specific coefficients 
(where k — 0, . . . , 6 and i = 1, . . . , f). The continuous red curve indicates the yield response 
as compared to the baseline scenario, i.e., the effect averaged over all crop models. The red dashed 
curves indicate the first and third quartiles of climate change effect computed over all crop models 
The blue line indicates zero yield difference. 


The fitted curves obtained for CO 2 response of wheat in Wageningen, the 
Netherlands are displayed in Figure 5b. Compared to maize, the effect of a CO 2 
concentration increase is stronger for wheat. This effect is highly variable among 
crop models; some models demonstrate positive slopes over the range of CO 2 con- 
centrations, whereas others show slopes close to zero. When averaged across aU of 
the crop models, the effect of die +2°C temperature increase is positive (yield gain) 
as soon as the CO 2 concentration increase reaches 100 ppm. 

Probabilistic thresholds 

An important advantage of the meta-model is that it handles the interpolation 
between temperature levels and between CO 2 concentration levels. Our meta-model 
can dins be used to calculate temperature and [CO 2 ] thresholds that lead to yield loss 
or yield gain. Table 1 shows the thresholds of [C02] increase required to obtain a 
pr(obability of maize and wheat yield gain higher than 0.5 (i.e., 50% chance of yield 
gain). These thresholds were computed for two values of temperature increase (-1-2 
and +4'’C) and four sites per crop. It would have been very difficult to compute them 
without a meta-model. Results show that, for wheat, climate change has 50% chance 
of resulting in a yield gain if [CO 2 ] increases by at least -1-117 ppm (depending on 
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Table 1. Levels of increased CO2 concentration (ppm) 
required to obtain a probability of maize and wheat yield 
gain higher than 0.5 (i.e., 50% chance of yield gain). Con- 
centrations were computed for two values of temperature 
increase (-1-2 and -|-4°C) and four sites per crop. Basehne 
scenario corresponds to [CO2] = 360 ppm. 



Temperature change 

Site 

-|-2“C 

-l-4''C 

Maize 

[CO 2 ] 

[CO 2 ] 

Rio Verde, Brazil 

>+360 

>+360 

Lusignan, France 

> -1-360 

> +360 

Morogoro, T^zania 

> +360 

> 4-360 

Ames, USA 

+269 

>+360 

Wheat 

[CO 2 ] 

[CO 2 ] 

Balcarce, Argentina 

+117 

+358 

Wongan Hills, Australia 

0 

+59 

New Delhi, India 

+112 

+278 

Wageningen, the Netherlands 

+83 

+222 


the site) and if temperature concurrently increases by +2°C (Table 1). Required 
levels of [CO 2 ] increases are much higher for maize (Table 1). Our meta-model has 
another advantage, although not fully illustrated here; it can be used to quantify the 
effects of temperature, [CO 2 ], and precipitation on yields, their interactions, and 
their variability between the four considered sites and between aU the considered 
crop models. It thus constitutes a powerful tool for exploring complex crop model 
responses to climate factors. 


Conclusions 

Our study shows that yield data simulated by an ensemble of complex dynamic- 
process crop models can be summarized by using statistical meta-models that 
are based on random coefficient regressions. These statistical models describe 
the between-crop model variability of the simulated yield data by using prob- 
ability distributions. They can be used to compute key quantities such as mean 
yield loss, percentiles of yield loss, and probabilities of yield loss as functions of 
temperature change and CO 2 concentration change. These statistical meta-models 
can thus be helpful for analyzing the risk of yield loss due to climate change at 
the locations where the original simulations were conducted. The capabilities of 
the meta-models were illustrated by using two yield datasets generated by the 
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AgMIP pilot model intercomparison studies for maize and wheat. The lesulto 
showed that: 


• ffigh CO 2 concentrations are likely to outweigh the negative effect of increasing 
temperature, leading to yield gains for wheat with climate change, but not for 
maize. 

• Thresholds of temperature increase and CO 2 concentrations that lead to a spe- 
cific risk level of maize and wheat yield loss can be computed by using the 
meta-model. 

. The divergence between file maize and wheat crop models and therefore the 
uncertainty in the simulated results increases as a function of temperature (the 
higher the temperature change, the higher the between-crop model variability). 

In the future, the proposed meta-model could be extended in two different ways. 
First, it could be applied to a dataset that includes simulations obtained for a hi gher 
number of sites. It will then be possible to include co-variables that describe site 
characteristics (e.g„ soil type, agricultural practices) in the meta-model in order to 
explain the origin of the between-site variability. Second, our meta-model could be 
extended in order to describe the between-year variability of yields and to analyze the 
risk of extreme yield values for different climate scenarios. This could be achieved 
by including one or several additional random effects in the meta-model and by 
fitting this model to yearly yield data. 
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